A New Method for Determining on the Lattice 
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We derive the two-point spectral correlation function of the Dirac operator with a specific external 
source in the e-regime of QCD. This correlation function has a unique and strong dependence on 
Fjy, and thus provides an novel way to extract Ft, from lattice simulations. We test the method in 
a quenched lattice simulation with staggered fermions. 
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INTRODUCTION 

One of the outstanding challenges to lattice gauge 
theory is the computation of physical observables that 
depend strongly on having very small quark masses in 
QCD. Here we exploit the fact that the low-lying spec- 
trum of the Dirac operator in finite volume is particularly 
sensitive to the observables of spontaneously broken chi- 
ral symmetry. It is well established U that the chiral con- 
densate E = can be extracted from measurements 
of the low-energy Dirac spectrum. We carry this program 
one step further to determine the next low-energy con- 
stant, the pseudoscalar decay constant i^^, with similar 
high precision. We show that a certain spectral correla- 
tion function of the Dirac operator depends on in a 
unique and quite spectacular way. Based on this depen- 
dence we propose and demonstrate a novel method for 
measuring F^^ in lattice gauge theory simulations. The 
method is general for systems with spontaneous break- 
ing of symmetries, and indeed the universal finite- volume 
scaling formulas have wide application in the context of 
condensed matter physics as well 

Conventionally, measurements of -FV on the lattice are 
carried out in the so-called p-regime using the 2-point 
function of the axial current. One then aims at lattices 
large enough that the Compton wavelength of the Gold- 
stone bosons is much smaller than the lattice size while 
still performing an extrapolation to the chiral limit of 
very light u and d quarks Q. The method we propose 
avoids such issues by going to the e-regime. 

The low-energy effective Lagrangian of the e-regime in 
QCD is dominated by the zero-momentum modes of 
the pseudo-Goldstone bosons. In the absence of external 
sources, the leading term in the associated e-expansion 
is proportional to the quark mass to. It has been sup- 
posed that needs to be computed either from the tiny 
;erturbative correction to the leading-order result for S 
Q or from an appropriate space-time correlation func- 
tion . Quenched Monte Carlo simulations have demon- 
strated that such a procedure is feasible, but numerically 



challenging j^]- 

Once external currents are included, the leading-order 
effective Lagrangian of the e-regime depends not only on 
E, but also on F^. We shall introduce an external vec- 
tor source that can be interpreted either as an imagi- 
nary chemical potential for isospin [lOl or as twisted 
boundary conditions for the gauge potentials 0. The 
advantage of an imaginary isospin chemical potential is 
twofold. First, the associated Dirac operator has a posi- 
tive definite determinant and thus becomes amenable to 
numerical simulations. In addition, with an imaginary 
isospin potential, the Dirac operator is anti-Hermitian 
and thus its eigenvalue spectrum lies entirely on the imag- 
inary axis. 

For non-zero baryon chemical potential, the -F^- 
dependent spectral correlation functions in the e-regime 
[12. 0, can be used to glean information about chi- 
ral symmetry breaking Tsj. However, the spectrum of 
the associated Dirac operator is complex, which makes 
the numerical determination of ^FV more demanding. 
Quenched results in that direction have recently been 
presented [lij . 

We consider instead the correlation function 

p{Xi,X2;iHiso) = ( ^<^(Ai - A„)^(5(A2 - A„) ^ 

- (5]5(Ai-A„)y^5(A2-A,„))(l) 



between the densities of eigenvalues iXn of the anti- 
hermitian operator where 

L'+l/'n = [0(A) + i^iisoTolV'n = ^AnVn, (2) 

and the eigenvalues iXm of the likewise anti-hermitian 
operator D_, defined by 



D^'lpn = [Ip[A) - i^lisolQ\i^n = iXm^n- 



(3) 



Here lp{A) is the Dirac operator associated with the 
gauge potential v4^. Considering the operators and 
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D- as acting on two separate flavors leads to a theory in 
which the w-quark has chemical potential +i^iso, while 
the d-quark has chemical potential —ifiiso- Using the ef- 
fective low energy theory for QCD in the e-regime we will 
calculate the correlation function in the microscopic 
limit where both mSV and fi^^^F^V are held fixed as 
the four-volume V is taken to infinity. We can measure 
E very accuratel y fr om the distributions of the smallest 
eigenvalues of Ip [ig • Then the spectral correlation func- 
tion m for eigenvalues of the order 1/{YjV) provides, as 
we shall see, a parameter-free determination of F^^ in the 
chiral limit. Computationally, all that is required is the 
determination of a modest number (typically, the first 
10-20 will suffice) of smallest eigenvalues of the Dirac op- 
erators D±. As an additional check, one can verify that 
the dependence is through the combination iJ-fg^F^V only. 

In this paper we will provide formulas and numeri- 
cal results relevant to quenched QCD. Shortcomings of 
quenching are well known, and we use this approxima- 
tion for illustrative purposes only. The case of dynamical 
quarks will be presented elsewhere. 



THE CORRELATION FUNCTION IN THE 
e-REGIME 

The first step is to derive the quenched susceptibility, 
defined as 

x{mi,m2;ifiiso)= lim \dmid,n2 log Z„(mi, m2; iAfiso)(4) 

n^O n 

where the limit n — > indicates use of the replica method 
(lif . In terms of the eigenvalues A and A the susceptibility 
is 



at fixed gauge field topology i/. The integration is over 
the Haar measure of U{2n), and we have defined 



X(mi,m2; j^iis 



^ iXn + "11 ^ fXrn + / 




The spectral correlation function ^ then follows from 
the discontinuity across the imaginary axis of both mi 
and 7Ti2, 

p(Ai, A2;i/iiso) = 7^Discx("ii,"i2;«Miso)|™i=iAi . (6) 

The QCD partition function Z„(mi, m2; i/iiso) entering 
Eq. 10} contains 2n fermions, half of which have mass 
nil and chemical potential i/iiso while the other half have 
mass 7712 and chemical potential — i/iiso- In the e-regime 
the leading term in the partition function is jll| 

Zn{■ml,m2;i^J,iso) ^ (7) 
jdU det(f7)''e^^^'''-°'^''[^"^l['^^^-^l+5^^'^'(*^*'^+*'^^*) 
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For n ~ 1 the partition function ((TJ reduces to 

Zi(mi,TO2;i^iso) (9) 

Jo 

and from this all partition functions for n > 2 can be 
obtained via 

(mim2)"("^^'Z„(TOi, m2; i^J.iso) (10) 
= Aidet [{■ml^rn^)''im2^m2yZliml,m2■,i^J,iso)] ■ 

Here Dn is a normalization factor and k,l = 0, 1, ..,n— 1. 

Recently it has been realized [l^Eoj that to obtain the 
correct replica limit rt — > in Eq. one can make use of 
the integrability relations satisfied by the partition func- 
tions. Equation H10|l has the structure of a r- function, 
implying that the Z„ satisfy the Toda lattice equation, 

Ati^V^Y.-^ ?7iit^mi m2dm2 log Z„ (mi , m2 ; i^iso) (11) 
= {mim2Y 



x2 Zn+i{mi, m2; i^iiso)Zn-ii'mi, m2; ifiiso) 



[Z„(toi, m2; i^iso)]^ 

where Dn fixes the coefficient on the Ihs 0| . Taking the 
n ^ limit of Eq. and comparing to Eq. (0} we find 



x{mi,m2;ifiiso) 
= 4mim2Zi(mi, 7712; i/iiso)-Z'-i(TOi, m2; i/iiso) 



(12) 



The ri — > limit in Eq. H12(l has naturally brought in 
the partition function with n = —1, i.e., one quark of 
bosonic statistics Like its fermionic analogue, the 

bosonic partition function is determined by the symme- 
tries of the underlying QCD Lagrangian. In the bosonic 
case, moreover, one must take care to ensure convergence 
of the partition function. As the purely imaginary chem- 
ical potential does not affect the hermiticity of the Dirac 
operator this does not lead to additional constraints (as 
opposed to the case of a real chemical potential |l3Lll4|'). 
The result for the bosonic partition fmiction is 

Z_i(?7ii,m2;«Aiiso) = e^^"^"^^=° (13) 

/OO 
dXXe-^^''^'''"'^^K^{XmiY.V)KA>^rn2T.V). 

The quenched susceptibility now follows from Eq. H12|l . 
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x(mi,TO2;z//iso) 

/oo 
dAAe"^^^'^'-^'if,.(ATOiSF)i4'^(ATO2S^). 

Taking the discontinuity as in Eq. (|SJ| we obtain the desired correlation function, 
p(a,|2;*Miso) = Ci6 /'dAAe2^^"''-^'j,(A6)J.(A|2) 



(14) 



(15) 



■ exp 



dX A e-'^^'^-^V.(Aei) J.(A6) 



where we have defined the scahng variables = XiT,V 
and ^2 = A2S1^. The spectacular change in this corre- 
lation function when /iiso is made non-zero can be seen 
in Fig. n A fit to Monte Carlo data with ^iso using 
Eq. (|15|l will then readily produce a measurement of Ft^ . 
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6 J.+i(6) J.(6) - 6^.+i(6)>/.(6 




FIG. 1: The correlation function with one eigenvalue fixed 
at ^2 = 4, for F^m^aVV = 0.0717 (full) and for /iiso = 
(dashed). The 5-function peak at = ^2 for /liso ~ has not 
been shown. 



The pronounced effect of /ijso 7^ in the region where A 
is close to A is not difficult to explain. For /iiso = both A 
and A are eigenvalues of Ip, and this leads to a (5- function 
contribution to the correlation function stemming from 
the non-compact integral in the bosonic partition func- 
tion (|14|l . which diverges for /Xiso = at Ai = A2. In this 
way we recover from Eq. (|15|l the known result |22| |. 



P(6,l2; JMiso = 0) 

= si^i - 6)|[j.'(6) - j.+i (6)^.-1 (a)] 



(16) 



The (5-function is not shown in the /iiso = curve in 
Fig. n When /iiso is non-zero, A is an eigenvalue of D+ 
while A is an eigenvalue of D- ; the effect of /iiso is there- 
fore to smooth out the 5-function into a pronounced peak 
for A near A. 



NUMERICAL SIMULATIONS 

To test the method, we have performed simulations 
of quenched QCD using staggered fermions for V — 8'^ 
(/iiso = 0.01) and V = 12^ (/iiso = 0.002). We have cho- 
sen to work with the standard Wilson plaquette action 
at P — b.7 and with conventional, unimproved Dirac op- 
erators. In this way we are sure to have no ambiguities in 
the identification of the coset space of spontaneous chiral 
symmetry breaking, which here is U(2n). The analy- 
sis presented above is based on the coset space SU{2n). 
We can account for the extra U{1) factor here by set- 
ting V — Q va. the formulas and comparing to numerical 
results without fixed topology (see for example the first 
paper of ref. 0). Simulation at weaker coupling or the 
use of improved actions and Dirac operators will induce a 
crossover from 1/ = behavior to an explicit dependence 
on topological index |25j . but it is not our purpose to 
explore that aspect here. 

We include a chemical potential on the lattice in the 
standard way [l^. For an imaginary chemical poten- 
tial this amounts to including a constant abelian gauge 
field with non-vanishing timelike component /x only. It 
leaves the hermiticity properties of the Dirac operator 
unchanged, allowing the computation of the lowest-lying 
eigenvalues with the Ritz variational algorithm [23| . 

When extracting physical observables we must keep 
in mind that the continuum theory describes 4 tastes of 
quark. In our simulation the staggered Dirac matrix is 4 
times larger than that of a single continuum quark. The 
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statistical properties of the eigenvalue spectrum of that 
matrix behave as in a theory of one species in a volume 
four times as big. Therefore, to determine the values 
of and E from the staggered eigenvalue spectrum we 
replace V in the analytical predictions by AV. 

First we measure E. We can do that by fitting individ- 
ual eigenvalue distributions to the analytical expressions 
|l7j . It is a non-trivial prediction that fc-point correlation 
functions of Z)+ and Z3_ separately are independent of 
/j,iso in the microscopic limit. This follows from Eq. ^ by 
taking B proportional to the unit matrix, and it ensures 
in particular that individual eigenvalue distributions are 
/.tiso-independent in this limit |2ll | . Alternatively, one can 
use the 2-point correlation functions of either D_|_ or D_ . 
As follows from the argument above, these 2-point func- 
tions are /iiso-independent, and they can provide indepen- 
dent determinations of E ||2^ . Using the latter approach 
for our 12^ data, a best fit gives in lattice units the bare 
value E = 0.0634(3) with x^/dof = 0.55. This is consis- 
tent with what we find by fitting either individual eigen- 
value distributions or the overall flat eigenvalue plateau. 
We can then determine by a fit of the measured cor- 
relation function to the analytical result (|15|l . In Fig. 
12 we show the data for V — 8* and compare them to the 
analytical curves using the best fit to from the 12'* 
data, as explained below. On V = 8'^ with P = 5.7 only 
the first few eigenvalues have distributions in agreement 
with the predictions of the e-regime. For this reason we 
base our measured values of E and on our 12^ lattice 
simulations. 



y = S' 10000 configs 
\^.J'y" = 0.159 



correlation function 



Pint (a;; ^^^iso) = 




FIG. 2: The correlation function with fixed ^2=4 mea- 
sured on our 8* lattice. The curves correspond to 11511 with 
F■^^Hso^/V = 0.159 (full) and /liso = (dashed). 



Requiring that the eigenvalue ^2 fall within one bin 
around a fixed value (here f 2 = 4) means that a large frac- 
tion of the lattice configurations is obviously not used. In 
order to improve the statistics we consider the integrated 



dCp(a; + ^,C;*Miso)- 



(17) 



In Fig. Owe show this integrated correlation function as 
measured on our 12^ ensemble. 

The best fit gives us — 0.1245(18) in bare lattice 
units with a /dof — 0.33. This value is consistent 
with the result = 0.118(7) of ref. [l^l, which uses the 
same action and gauge coupling /3 = 5.7 (but of course a 
different method for extracting F^r). 

It is well known that the unimproved staggered action 
leads to serious scaling violations. To quote a result in 
physical units would, at this value of /3, require the use 
of improvement methods to reduce lattice artifacts. This 
is beyond the scope of this paper. 
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FIG. 3: The integrated correlation function l)17|l with ^min = 
4 and ^max ~ 43 . The value F^/iisoV'V' = 0.0717 is the best 
fit. 



SUMMARY 

We have proposed and tested a new method to mea- 
sure Ftt in lattice gauge theory simulations, using the 
pronounced F,r-dependence of a specific 2-point correla- 
tion function of the Dirac operator in the e-regime. 

We have performed quenched lattice simulations that 
demonstrate excellent agreement with our analytical pre- 
dictions. Our study illustrates the ease with which this 
method can be implemented, and the high precision that 
can be achieved. 
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